setwd("E:\\5hmc_file\\2_5hmc_yjp_bam\\ASM\\bayes_pvalue_beta0")
group1=c("X2B_X1T","M8_M7","M6_M5","M2_M1","M48_M47","M50_M49","M28_M27","M30_M29","M26_M25","M35_M36","M18_M17","M20_M19","M22_M21","M40_M39")

i=1
fn1=paste0(group1[i],".bayes_p.txt")
file=read.table(fn1,head=T,sep = "\t")
file$unitID=paste(file$chrom,file$position,file$ref,file$var,sep=":")
file=file[file$normal_bayes_pvalue<0.05|file$tumor_bayes_pvalue<0.05,]
file=data.frame(unitID=file$unitID,normal_bayes_beta0=file$normal_bayes_beta0,normal_bayes_pvalue=file$normal_bayes_pvalue,tumor_bayes_beta0=file$tumor_bayes_beta0,tumor_bayes_pvalue=file$tumor_bayes_pvalue)
file$normal_group="nosig"
file[file$normal_bayes_beta0>0&file$normal_bayes_pvalue<0.05,]$normal_group="up"
file[file$normal_bayes_beta0<0&file$normal_bayes_pvalue<0.05,]$normal_group="down"
file$tumor_group="nosig"
file[file$tumor_bayes_beta0>0&file$tumor_bayes_pvalue<0.05,]$tumor_group="up"
file[file$tumor_bayes_beta0<0&file$tumor_bayes_pvalue<0.05,]$tumor_group="down"
file$pattern=paste0("normal_",file$normal_group,"-tumor_",file$tumor_group)
names(file)=c("unitID",paste(group1[i],names(file)[2:5],sep="_"),names(file)[6:8])

normal_up_tumor_nosig=file[file$pattern=="normal_up-tumor_nosig",][,1:5]
normal_up_tumor_up=file[file$pattern=="normal_up-tumor_up",][,1:5]
normal_down_tumor_nosig=file[file$pattern=="normal_down-tumor_nosig",][,1:5]
normal_nosig_tumor_up=file[file$pattern=="normal_nosig-tumor_up",][,1:5]
normal_down_tumor_up=file[file$pattern=="normal_down-tumor_up",][,1:5]
normal_up_tumor_down=file[file$pattern=="normal_up-tumor_down",][,1:5]
normal_down_tumor_down=file[file$pattern=="normal_down-tumor_down",][,1:5]
normal_nosig_tumor_down=file[file$pattern=="normal_nosig-tumor_down",][,1:5]
###单发样本取交集，upup有14个交集，downdown8个，其余的没有交集。
for (i in 2:6) {
fn1=paste0(group1[i],".bayes_p.txt")
file=read.table(fn1,head=T,sep = "\t")
file$unitID=paste(file$chrom,file$position,file$ref,file$var,sep=":")
file=file[file$normal_bayes_pvalue<0.05|file$tumor_bayes_pvalue<0.05,]
file=data.frame(unitID=file$unitID,normal_bayes_beta0=file$normal_bayes_beta0,normal_bayes_pvalue=file$normal_bayes_pvalue,tumor_bayes_beta0=file$tumor_bayes_beta0,tumor_bayes_pvalue=file$tumor_bayes_pvalue)
file$normal_group="nosig"
file[file$normal_bayes_beta0>0&file$normal_bayes_pvalue<0.05,]$normal_group="up"
file[file$normal_bayes_beta0<0&file$normal_bayes_pvalue<0.05,]$normal_group="down"
file$tumor_group="nosig"
file[file$tumor_bayes_beta0>0&file$tumor_bayes_pvalue<0.05,]$tumor_group="up"
file[file$tumor_bayes_beta0<0&file$tumor_bayes_pvalue<0.05,]$tumor_group="down"
file$pattern=paste0("normal_",file$normal_group,"-tumor_",file$tumor_group)
names(file)=c("unitID",paste(group1[i],names(file)[2:5],sep="_"),names(file)[6:8])

normal_up_tumor_nosig=merge(normal_up_tumor_nosig,file[file$pattern=="normal_up-tumor_nosig",][,1:5],by="unitID")
  normal_up_tumor_up=merge(normal_up_tumor_up,file[file$pattern=="normal_up-tumor_up",][,1:5],by="unitID")
  normal_down_tumor_nosig=merge(normal_down_tumor_nosig,file[file$pattern=="normal_down-tumor_nosig",][,1:5],by="unitID")
  normal_nosig_tumor_up=merge(normal_nosig_tumor_up,file[file$pattern=="normal_nosig-tumor_up",][,1:5],by="unitID")
  normal_down_tumor_up=merge(normal_down_tumor_up,file[file$pattern=="normal_down-tumor_up",][,1:5],by="unitID")
  normal_up_tumor_down=merge(normal_up_tumor_down,file[file$pattern=="normal_up-tumor_down",][,1:5],by="unitID")
  normal_down_tumor_down=merge(normal_down_tumor_down,file[file$pattern=="normal_down-tumor_down",][,1:5],by="unitID")
  normal_nosig_tumor_down=merge(normal_nosig_tumor_down,file[file$pattern=="normal_nosig-tumor_down",][,1:5],by="unitID")
}
结果在20201104下的intersect后缀文件中。
write.table(normal_up_tumor_up,"../20201104_at.least.one.is.AShM.beta0.pic/normal_up_tumor_up.intersect.txt",quote = F,row.names = F,sep="\t")
write.table(normal_down_tumor_down,"../20201104_at.least.one.is.AShM.beta0.pic/normal_down_tumor_down.intersect.txt",quote = F,row.names = F,sep="\t")

										###22个位点的beta0可视化
library(dplyr)
library(ggplot2)
normal_up_tumor_up$group="normal_up_tumor_up"
normal_down_tumor_down$group="normal_down_tumor_down"
file=rbind(normal_up_tumor_up,normal_down_tumor_down)

file$mean.normal.beta0=rowMeans(file[,seq(2,24,4)])
file$mean.tumor.beta0=rowMeans(file[,seq(4,24,4)])
std <- function(x) sd(x)/sqrt(length(x))###算标准误差
file$std.normal.beta0=apply(file[,seq(2,24,4)], 1, std)
file$std.tumor.beta0=apply(file[,seq(4,24,4)],1,std)
file$normal.beta0.upper=file$mean.normal.beta0+file$std.normal.beta0
file$nomal.beta0.lower=file$mean.normal.beta0-file$std.normal.beta0
file$tumor.beta0.upper=file$mean.tumor.beta0+file$std.tumor.beta0
file$tumor.beta0.lower=file$mean.tumor.beta0-file$std.tumor.beta0
anno=read.csv("../20201104_at least one individuals is AShM/intersect.DC.anno.hg19_multianno.csv",head=T)
anno$unitID=paste(anno$Chr,anno$Start,anno$Ref,anno$Alt,sep=":")
file=merge(anno,file,by="unitID")

normaldata=select(file,unitID,mean.normal.beta0,normal.beta0.upper,nomal.beta0.lower,avsnp150,group)
normaldata=arrange(normaldata,mean.normal.beta0,group)
normaldata$num=1:dim(normaldata)[1]
orders=select(normaldata,unitID,num)

p1=ggplot(normaldata,aes(y=mean.normal.beta0,x=num))+scale_x_continuous(breaks = normaldata$num,labels = normaldata$avsnp150,position="top")+
  geom_point(aes(color=group),size=1.5)+geom_errorbar(aes(ymin=nomal.beta0.lower,ymax=normal.beta0.upper,color=group),width=0.5)+labs(x="",y="con β0")+
  geom_hline(yintercept = 0,color="red",linetype="dashed")+scale_color_manual(values = c('#708090','red'))+coord_flip()+theme_classic(base_size = 15)+
  theme(panel.grid.minor = element_blank(),panel.border = element_blank(),axis.line =element_line(size = 0.5))+theme(legend.position="none")

tumordata=select(file,unitID,mean.tumor.beta0,tumor.beta0.upper,tumor.beta0.lower,avsnp150,group)
tumordata$na=""
tumordata=merge(tumordata,orders,by="unitID")

p2=ggplot(tumordata,aes(y=mean.tumor.beta0,x=num))+scale_x_continuous(breaks = tumordata$num,labels = tumordata$na)+
  geom_point(aes(color=group),size=1.5)+geom_errorbar(aes(ymin=tumor.beta0.lower,ymax=tumor.beta0.upper,color=group),width=0.5)+
  labs(x="",y="case β0")+geom_hline(yintercept = 0,color="red",linetype="dashed")+scale_color_manual(values = c('#708090','red'))+theme_classic(base_size = 15)+
  theme(panel.grid.minor = element_blank(),panel.border = element_blank(),axis.line = element_line(size = 0.5))+coord_flip()+theme(legend.position="none")

layout <- matrix(c(1, 1, 1, 1, 2, 2, 2), nrow = 1)
multiplot(plotlist=list(p1,p2),layout=layout)